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The  authors  would  like  to  thar.k  C,  E.  Harbec*  o i  the  U.  r' ,  Geologi* 
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SUMMARY 


The  problem  of  predicting  downstream  changes  4n  river  temperatures 
due  to  the  addition  of  various  forms  of  thermal  pollution  can  be  expected 
to  assume  greater  importance  as  man  increasingly  uses  river  water  as  a 
coolant  and  *s  a  medium  for  transporting  biological  A  possible 

result  of  such  pollution  is  the  maintenance  of  ice-free  river  reaches 
during  the  winter  which  might  have  beneficial  consequences  to  river  traf¬ 
fic  on  navigable  waterways.  Previous  methods  for  predicting  downstream 
temperature  decrease,  which  could  be  used  to  estimate  lengths  of  ice- 
frre  reaches,  are  limited  in  that  only  a  few  of  the  Several  pertinent 
meteorological  factors  determining  heat  loss  rates  are  considered. 

This  report  formulates  a  differential  equation  far  the  steady-state 
heat  balance  of  a  volume  clement  of  a  river,  which  leads  to  the  expres¬ 
sion 

Twx 

r  d  T  w 
x  <-x  y  q* 

T  wo 

where  x  is  the  distance  downstream  from  the  thermal  pollution  site  to 
the  cross  section  where  water  temperature  equals  T^,  Tvv#  is  the  water 
temperature  at  x  equals  zero,  Q+  is  the  rate  of  heat  loss  from  the  water 
surface,  and  Cx  is  a  constant  which  includes  flow  velocity  and  depth, 

Q*  is  the  sum  of  heat  losses  duo  to  evaporation,  convection,  long-  and 
short-wave  radiation,  and  other  processes,  each  of  which  is  evaluated  by 
an  empirical  or  theoretical  expression.  For  given  hydrometric  and 
meteorological  conditions  (air  temperature,  wind  speed,  humidity,  cloud 
cover,  cloud  height,  etc),  I  /Q*  is  a  complicated  function  of  water  tem¬ 
perature.  A  computer  program  numerically  evaluates  the  integral  and 
calculates  the  distance  x  for  given  values  of  Twx<  The  value  of  x  when 
Twx  equals  OC  is  taken  as  the  length  of  the  ice-free  reach. 

There  are  two  principal  limitations  to  accurate  prediction  of  down¬ 
stream  temperate  re  changes  .  Fi  rst,  it  is  difficult  to  estimate  the  de¬ 
gree  of  lateral  mixing  in  a  natural  river.  Second,  there  is  considerable 
uncertainty  as  i  o  how  t:>  best  o valuai  e  convective  and  evaporative  heat 
losses  under  conditions  of  local  atmospheric  instability.  Twa  approaches 
are  used  in  this  report,  both  based  on  the  application  of  (he  Bowen  ratio 
to:  (I)  the  widely  used  Kohler  equation  for  evaporative  heat  loss;  and 
(4)  an  equation  for  convective  heat  loss  developed  empirically  by  Russian 
workers  under  wintertime  conditions.  The  latter  equations  predict  quite 
well  the  ebserved  lengths  of  ice-free  reaches  associated  with  power  and 
sewage  plants  along  the  Mississippi  River  at  Minneapolis  -  St.  Paul, 

M  mnesota. 

It  is  shown  that  significant  portions  of  the  St,  Lawrence  Seaway  can 
be  kept  ice-free  by  the  installation  of  nuclear  reactors  at  appropriate 
locations. 


THE  EFFECTS  OF  THERMAL  POLLUTION  ON  RIVER  ICE  CONDITIONS 
I.  A  GENERAL  METHOD  OF  CALCULATION 

by 

S.  L.  Dingman,  W.F.  Weeks,  and  Y.C.  Yen 
INTRODUCTION 

The  general  problem  of  the  magnitude  and  extent  of  downstream  changes 
in  river  temperature  due  to  the  addition  of  "thermal  pollution"  in  various 
forms  is  assuming  greater  importance  as  man  increasingly  uses  river  water 
as  a  coolant  and  as  a  medium  for  transporting  biological  wastes.  This  paper 
treats  a  special  case:  the  effects  of  such  heat  additions  on  river  ice  conditions. 
The  basic  approach,  however,  is  applicable  to  the  general  problem. 

If  temperature  increases  due  to  thermal  pollution  are  of  sufficient  magni¬ 
tude,  an  appreciable  reduction  in  thickness,  or  even  complete  disappearance, 
of  the  winter  ice  cover  will  result.  The  thermonuclear  production  of  electrical 
power  with  the  associated  production  of  large  quantities  of  waste  heat  is  now 
becoming  competitive  with  the  more  conventional  methods  of  power  production. 
If  it  can  be  established  that  this  waste  heat  is  sufficient  to  keep  major  portions 
of  navigable  rivers  free  from  ice,  the  proper  positioning  of  the  thermonuclear 
reactors  might  permit  the  effective  use  of  the  river  by  shipping  throughout  a 
significant  portion  of  the  winter.  The  advantages  of  being  able  to  calculate  the 
effects  of  changes  in  the  amount  of  thermal  pollution,  the  hydrometric  param¬ 
eters  of  the  stream,  and  the  weather  conditions  on  the  position  of  the  ice  front 
below  a  pollution  site  are  therefore  obvious. 


PREVIOUS  WORK 

While  heat  balance  studies  of  lakes  are  numerous,  very  few  such  investi¬ 
gations  exist  for  rivers.  Eckel  and  Reuter  (IV50),  who  recognized  the  need 
for  such  work,  were  concerned  primarily  with  predicting  temperatures  in 
rivers  under  natural  conditions.  Therefore  they  attempted  to  develop  an  equa¬ 
tion  to  estimate  the  changes  in  temperature  of  a  reach  of  river  or  moving 
parcel  of  water  with  time,  rather  than  with  distance.  They  considered  heat 
budget  items  associated  with  short-  and  long-wave  radiation  and  evaporation 
and  convection  to  the  atmosphere.  Because  of  the  complexity  of  the  total  heat 
balance  term,  the  differential  equation  they  developed  could  not  be  solved  in 
closed  form,  and  they  used  graphical  iterative  procedures  which  gave  satis¬ 
factory  agreement  between  measured  and  calculated  river  temperatures,  at 
least  for  clear  summer  days. 

Pruden  et  al.  (1954)  studied  the  heat  balance  of  the-  St.  Lawrence  River 
in  an  attempt  to  determine  the  feasibility  of  keeping  much  of  the  river  ice- 
free  for  parts  of  the  winter.  They  considered  the  meteorological  dependence 
of  the  various  heat  balance  terms  in  somewhat  more  detail  than  had  Eckel  and 
Reuter.  Because  plots  of  the  heat  loss  rate  (Q^)  vs  air  temperature  (Ta)  and 
vs  water  temperature  (Tw)  were  close  to  linear  (using  average  monthly  values 
for  meteorologic  conditions),  they  wrote? 

Q*  =  183.6  +  15. 5T  +  43.1  (T  -  T  ) 

a  w  a 


(1) 
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where  Q*  is  in  cal  cm*2  day"1  and  T^  and  T,  are  in  degrees  C.  This  ex¬ 
pression  for  Q*  was  then  incorporated  in  a  differential  equation  for  cooling 
rate  with  distance,  which  is  similar  to  the  one  derived  in  this  paper.  The 
inverse  of  their  expression  for  Q*  can  be  integrated  in  closed  form,  and 
they  found  a  simple  exponential  decrease  of  water  temperature  with  distance 
which  is  very  similar  to  that  found  in  the  present  paper  (Figure  3  shows  a 
typical,  nearly  linear  exponential  curve).  Their  approach  is,  however,  less 
satisfactory,  as  their  expression  for  Q#  was  developed  for  monthly  average 
values  at  one  general  location  only,  and  does  not  take  into  account  the  many 
meteorologic  variables  which  affect  the  heat  loss  rate.  Wind  speed,  for 
example,  is  shown  here  to  be  an  extremely  important  factor  in  determining 
Q*. 

Ince  and  Ashe  (1964)  used  observations  from  the  St.  Lawrence  River  to 
compare  two  methods  of  calculating  downstream  water  temperature  de¬ 
creases.  Both  of  these  arc  based  on  a  finite  difference  approach  for  a  parcel 
of  water  moving  at  the  average  velocity  of  the  river: 


T  T  -  Qt 

W2  Twl  "  DPwCp 


(2) 


where  Qt  is  heat-loss  rate  (cal  sec"1 ),  D  is  river  discharge  (cm1  sec-1 ),  p 
is  mass  density  of  water  (g  cm*1),  Cp  is  specific  heat  of  water  (calg*1  'C*1 ), 
and  (Tw,  *  Tv/[ )  is  temperature  change  of  the  water  parcel  over  some  time 
period  (^C).  In  the  first  method,  heat  exchanges  due  to  radiation,  evapora¬ 
tion,  convection  and  precipitation  were  accounted  for  in  evaluating  Q*,  using 
daily  averages  of  the  meteorological  variables,  but  no  details  of  the  various 
expressions  were  presented.  This  approach  was  compared  by  Ince  and  Ashe 
(1964)  to  that  using  the  simple  formula 


Q*  =  46.4  (Tw  -  Ta)  (3) 

recommended  by  the  Canadian  Joint  Board  of  Engineers,  where  Q*  is  heat 
loss  rate  in  cal  cm'2  day"1  and  (Tw  -  Ta)  is  based  on  daily  averages  of  air 
and  water  temperatures  in  °C.  Both  approaches  gave  satisfactory  estimates 
of  water  temperature  under  certain  conditions. 

While  the  basic  approaches  of  these  earlier  studies  are  sound,  and  led 
to  some  success  in  river  temperature  prediction,  they  are  all  limited  in 
some  respect.  The  principal  limitations  lie  in  the  formation  of  the  expres¬ 
sions  for  heat  loss  rate  Q*.  The  present  paper  largely  overcomes  these 
limitations  by  use  of  formulas  for  the  heat  budget  components  that  can  be 
used  with  daily  averages  of  meteorological  variables  and  by  use  of  a  com¬ 
puter  program  which  allows  numerical  integration  of  a  more  complete  and 
complicated  expression  for  1/Q*.  The  limitations  of  the  present  approach 
arise  from:  1)  uncertainties  in  estimating  the  evaporative  and  convective 
components  of  heat  loss  from  open  water  during  periods  of  pronounced  local 
atmospheric  instability;  2)  difficulties  in  evaluating  the  degree  of  mixing  of 
the  heated  water  through  the  stream  cross-section;  and  3)  the  assumption 
of  steady  conditions. 
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usr  of  symbols 

When  units  used  in  the  derivation  in  the  text  are  different  Ifom  those  used 
m  t tie  e computer  program,  t)u-  units  used  in  the  computer  program  are  under¬ 
lined, 


Text 

Meaning 

Units 

Fortran 

A 

Snow  accumulation  rate 

g  cm **  day*1 

ACC 

C 

Cloud  cover 

dimensionless 

clix: 

Ci 

Specific  heat  of  ice 

cal  g->  *C*‘ 

C 

P 

Specific  heat  of  water 

cal  g*'  ’C*' 

CP 

D 

HU  er  discharge 

cm1  sec"1 

E 

Evaporation  rate 

cm  day'1 

H 

Cloud  height 

m 

CL  Dll 

J 

Thermal  equivalent 

dy  cm  cal'1 

L 

Latent  heat  of  fusion  uf  ice 

cai  g'1 

P 

Atmospheric  pressure 

mb 

DUES,  P 

Qa 

Heat  added  to  river  by 
long-wave  radiation 

cal  cm‘!  day'1 

QA 

Qar 

Long-wave  radiation  reflected 
at  river  surface 

cat  cm'1  day'1 

OAR 

°B 

Net  heat  lost  from  river  as 
long-wave  radiation 

cal  cm  '4  day 

OB 

Qb« 

Long -wave,  back  radiation 
from  river  surface 

cai  cm'*  day'1 

0 13  S 

QCL 

Incoming  clear  sky  short-wave 
radiation 

cal  cm’!  day*1 

QCL 

qe 

Heat  lost  from  river  by 
evaporation 

cai  cm  ~l  day ' 1 

QE 

qf 

Heat  produced  by  fluid  friction 

cai  cm'4  day'1 

qg 

Heat  added  to  river  by  geo¬ 
thermal  heat  flow 

cal  cm'4  day*1 

qgw 

Heat  added  to  river  by  ground- 
water  inflow 

cal  cm*4  day  * 1 

qh 

Heat  lost  from  river  by 
convection  to  atmosphere 

cal  cm  "4  day*1 

OH 

Qp 

Heat  added  to  river  as  thermal 
pollution 

cal  sec *l 

POLL 

qr 

Heat  added  to  river  by  short¬ 
wave  radiation 

cal  cm'4  H.ty  *■ 

QH 
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LIST  OF  SYMBOLS  (Confd) 


Text 

Meaning 

Units 

Fortran 

°RI 

Short-wave  radiation  incident 
at  river  surface 

cal  cm*4  day*1 

QRI 

qrr 

Short-wave  radiation  reflected 
at  river  surface 

cal  cm*4  day** 

QRR 

°s 

Heat  lost  from  river  due  to 
in-falling  snow 

cal  cm'4  day*’ 

OS 

o, 

Heat  advected  into  elemental 
volume 

cal  sec*1 

Heat  conducted  into  elemental 
volume 

cal  sec"1 

Q» 

Heat  advected  from  elemental 
volume 

cal  sec*1 

a< 

Heat  conducted  from  elemental 
volume 

cal  sec*1 

Ql 

Heat -loss  rate  from  water 
surface 

cal  sec*1 

Q# 

Heat -  loss  rate  per  unit  area 
of  water  surface 

caltm*’  sec*1 
cal  cm*4  day*1 

QSTAR 

R 

Bowen  ratio 

dimensionless 

S 

Slope  of  river  surface 

dimensionless 

Ta 

Air  temperature 

*k,  *_£ 

ta,  TAC 

Ti 

Initial  water  temperature  pro¬ 
duced  by  thermal  pollution 

’K,  '_C 

TWIN,  TWINC 

T 

nat 

Natural  water  temperature 
abt  ve  thermal  pollution  site 

•K,  "C 

TWNAT,  TWNATC 

To 

Datum  temperature 

*  K _ 

Tw 

Water  temperature 

*K,  *_C 

TW,  TWC 

Tx 

Water  temperature  at  distance 
x  below  thermal  pollution  site 

*K,  *_C 

TWINIT,  TWFINL 

U 

Relative  humidity 

percent 

U 

V 

Visibility 

km 

V 

X 

Total  length  of  ice -free  reach 

cm,  m 

a 

Empirical  constant  (eq  20) 

A 

b 

Empirical  constant  (eq  20) 

B 

ea 

Vapor  pressure  of  air 

mb 

EA 

c  sa 

Saturation  vapor  pressure  of 
air  at  temperature  Ta 

mb 

VAP(TA) 
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LIST  OF  SYMBOLS  (Cont'd) 


Text 

Meaning 

Units 

Fortran 

e 

sw 

Saturation  vapor  pressure  of 
air  at  temperature  Tw 

mb 

V.AP(TW) 

Average  depth  of  river 

cm,  _m 

H 

k 

Thermal  conductivity  of  water 

ealem*1  sec"1  *C‘* 

V 

k 

Wind  velocity 

m  sec*1 

VA 

V 

w 

Average  river  flow  velocity 

cm  sec*1 ,  m  sec*1 

VW 

w 

Average  river  width 

cm,  m 

w 

/■ 

Distance  downstream  from 
thermal  pollution  site 

cm,  m 

X 

V 

Weight  density  of  water 

dy  cm*1 

y 

Ratio 

dimensionless 

FUG,  K 

Eddy  diffusivity  of  heat  in  air 

cm*  sec*1 

Kw 

Eddy  diffusivity  of  water  vapor 
in  air 

cm1  sec*1 

Pw 

Mass  density  of  water 

g  cm*1 

<r 

Stefan -Boltzmann  constant 

cal  cm*1  day*1*  K"4 

SBC 

BASIC  EQUATION 

The  differential  equation  fo^  i  vnstream  changes  in  water  temperature 
under  certain  simplifying  assun.pt  ms  is  developed  as  follows.  Referring 

to  Figure  1,  consider  a  stream,  carrying  a  steady  flow,  whose  geometry  La _ 

^eonstahfTrT space  and  time.  Assume  further  that  meteorological  conditions 
above  the  stream  are  spatially  and  temporally  constant.  If  water  tempera¬ 
ture  changes  in  the  downstream  (x)  direction  only,  the  amount  of  heat  ad- 
vected  into  volume  element  w-h-dx  in  unit  time,  Qj  (cal sec'1 ),  is 

Q»  °  PwvwhwVTw  -  T nat>-  (4) 

The  heat  transfer  to  the  volume  by  conduction,  Qt  ,  is 
dT 

Qa  =  -khw  .  (5) 

Here  pw  is  the  mass  density  of  the  water  (g  cm*1),  vw  is  the  flow  velocity 
of  the  water  (cm  sec"1 ),  h  is  the  flow  depth  (cm),  w  is  the  stream  width  (cm), 
Cp  is  the  specific  heat  of  the  water  (cal  g"1  *K'1),  considered  constant  in  the 
temperature  range  Tw  to  Tnat,  Twis  the  wafer  temperature  (*K),  k  is  the  eddy 
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Figure  1.  Definition  sketch  for 
basic  equation. 


thermal  conductivity  of  the  water  (cal 
cm*1  sec-l  'K*1)*  and  Tnat  is  the  natur' 
al  water  temperature  just  upstream 
from  the  thermal  pollution  site  (*K). 
The  amount  of  heat  advected  from  this 
volume  is 


(6) 

and  that  leaving  by  conduction  is 

rdTv,  d*  Tw  i 

q4  *  -khwLjk-  +  -g^-dxj. 

(7) 


The  water  surface  acts  as  a  heat  sink,  and  the  total  heat  lost  from  it.  fit  > 

is 


Ql  «  Q*-  w-  dx, 


(B) 


where  Q*  is  the  heat  loss  rate  per  cm2  of  water  surface.  A  heat  balance 
statement  for  the  volume  of  water  may  be  written  as 


Qj  +  Q2  =  Qj  +  Q4  +  Qt . 

which  reduces  to 


(9) 


dT 


Q*  =  -p  v  hC  — i— -  +  kh  j-3 -  . 

Kw  w  p  dx  dx* 

Neglecting  the  second  order  (conduction)  term,  eq  10  simplifies  to 


d2  T 


w 


(10) 


x  *  -p  v  hC 
rw  w  p 


dTw 

Tpr> 


di) 


where 


Ti 


T  -f  ^ 
nat  p  v  hwC 
rw  w  p 


(12) 


Tx  *  water  temperature  at  x. 


Tnat  28  "atural  water  temperature  just  upstream  from  the  thermal 
pollution  site,  and  Qp  is  the  heat  of  pollution  (cal  sec*1).  By  evaluating 
eq  11  at  different  values  of  Tx,  the  cooling  curve  downstream  from  the 
thermal  pollution  site  can  be  obtained. 
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Before  examining  Q*  further,  the  assumptions  involved  in  the  derive* 
lion  of  tq  11  should  be  considered.  In  most  cases,  the  assumption  of  steady 
flow  is  reasonably  good  for  ice*covered  rivers,  and,  unless  large  tributaries 
enter,  discharge  can  be  considered  nearly  constant  over  fairly  long  reaches. 
Where  large  tributaries  are  present,  or  where  there  are  significant  changes 
in  stream  geometry  along  the  reach  of  interest,  the  total  reach  can  be  sub* 
divided  and  each  subreach  considered  separately.  A  similar  approach  could 
be  used  to  consider  spatially  varying  meteorologic  conditions. 

The  assumption  that  water  temperature  is  constant  throughout  any  cross 
section  of  the  river  is  more  difficult  to  evaluate.  The  intensity  of  turbulent 
mixing  in  a  given  direction  equals  the  product  of  the  coefficient  of  turbulent 
diffusion  and  the  velocity  gradient  in  that  direction.  Thus  vertical  mixing 
tends  to  be  more  complete  (and  hence  temperatures  more  uniform  vertically) 
in  a  shallow,  fast  stream  than  in  a  deep,  slow  one.  By  the  same  reasoning, 
one  would  expect  narrow,  fast  streams  to  have  more  complete  lateral  mixing 
and  more  uniform  temperature  horisontally  than  wide,  slow  streams.  Sayre 
and  Chamberlain  (1964)  found  that  the  coefficient  of  lateral  turbulent  diffu¬ 
sion  is  approximately  proportional  to  the  3/2  power  of  the  mean  flow  depth, 
so  that  one  may  add  that,  velocity  gradients  and  other  things  equal,  greater 
depth  of  flow  promotes  lateral  mixing. 

If  there  is  a  basis  for  estimating  the  lateral  temperature  gradient  at 
various  points  along  the  reach  of  interest,  consecutive  sub  ’aches  of  vary¬ 
ing  width  can  be  used  to  approximate  the  natural  conditions.  Where  no 
such  information  exists,  one  must  be  satisfied  with  the  assumption  of  com¬ 
plete  isothermality  in  all  cross  sections.  For  an  example  of  how  erroneous 
this  assumption  can  be,  see  Mackay  (1966), 


EVALUATION  OF  Q* 

To  solve  eq  11,  the  heat  transferred  from  the  water  to  the  air  must  be 
written  as  a  function  of  water  temperature  Tw.  To  do  this,  a  second  heat 
balance  equation  is  written  in  which  all  the  modes  by  which  heat  can  enter 
and  leave  the  volume  of  water  w.h«dx,  aside  from  advection,  are  enumer¬ 
ated.  The  magnitude  of  the  heat  transferred  by  each  mode  is  written  as  a 
function  of  hydrometric  and  meteorological  variables,  and  the  sum  of  these 
equated  to  Q*.  As  shown  below,  the  terms  which  contribute  significantly  to 
Q*  can  be  evaluated  under  given  climatic  and  flow  conditions. - 

Q*  can  thus  be  written  as 

Q*  ■  Qr  *  qb  -  qe  *  qh  '  Qs  +  Qc  +  qgw  +  qf  (l  3) 

where 

QR  is  heat  added  at  the  river  surface  as  short-wave  radiation; 

Qg  is  net  loss  of  heat  by  the  exchange  of  long -wave  radiation  with 

the  atmosphere; 

Q^.  is  heat  loss  due  to  evaporation; 
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Q|.  is  sensible  heat  lost  by  conduction  and  subsequent  convection 
from  the  water  to  the  atmoaphere; 

Qg  ia  heat  lost  as  a  result  of  anow  falling  into  the  water; 

Qr  ia  heat  added  by  the  flow  of  geothermal  heat  through  the  stream 
°  bed; 

Qrw  is  heat  added  by  the  flow  of  ground  water  into  the  stream  through 
its  bed  and  banka; 

Q_  ia  heat  increase  due  to  fluid  friction. 

F 

With  these  definitions.  Q*  is  negative  when  there  is  a  net  loss  from  the 
volume  w>  h*  dx. 

°R 

Direct  measurements  of  solar  and  atmospheric  short-wave  radiation  at 
the  site  under  consideration  are  obviously  best  for  determining  the  magnitude 
of  the  energy  supplied  by  this  means.  Where  these  do  not  exist,  measure¬ 
ments  from  stations  at  some  distance  from  the  site  can  be  substituted,  pro¬ 
vided  there  is  reason  to  believe  that  the  cloud  cover  conditions  at  the  two 
sites  are  comparable  during  the  period  of  concern.  If  direct  measurements 
of  short-wave  radiation  are  lacking,  estimates  of  its  magnitude  can  be  made 
by  using  tables,  graphs,  or  formulas  relating  percentage  of  possible  sun¬ 
shine  or  cloud  cover  to  the  percent  of  possible  solar  radiation  which  reaches 
the  ground  (see,  for  example,  List,  1963,  Tables  151,  1 52;  Anderson,  1954, 
p.  74-77;  Jensen  and  Haiae,  1963,  App.  I).  Bolsenga  (1964)  has  provided 
tables  of  daily  sums  of  cloudless  sky  short-wave  radiation  for  various  lati¬ 
tudes  and  air  mass  conditions  which  can  be  used  as  a  basis  for  such  esti¬ 
mates.  For  the  purposes  of  the  present  paper,  the  amount  of  incident  short¬ 
wave  radiation  Qqt  is  estimated  using  the  relation  initially  suggested  by 
Angstrom  (List,  T963): 

°RI  3  °CL  l°*  35  +  0-061  (1°  -  c>]  (14) 

where  Q^l  *8  incoming  short-wave  radiation  for  a  cloudless  sky  and  C 
is  the  cloudiness  expressed  in  tenths  (i.  e, ,  with  complete  cloud  cover,  C  = 

A  portion  of  the  short-wave  radiation  which  reaches  the  water  surface 
is  reflected,  and  this  must  be  subtracted  from  the  incident  radiation  to  find 
the  energy  actually  absorbed.  Studies  at  Lake  Hefner,  Oklahoma  (Anderson, 
1954,  p.  78-88)  showed  that  water  surface  albedo  ranged  from  0,05  to  0,  30, 
was  a  strong  function  of  solar  altitude,  and  depended  to  some  extent  on  cloud 
conditions.  Unfortunately,  it  iB  difficult  to  use  Anderson's  relations  because 
the  empirical  constants  must  be  changed  to  correspond  with  the  type  of  cloud 
cover.  Koberg  (1964)  has  suggested  a  way  around  this  problem  by  plotting 
curves  for  reflected  short-wave  radiation  versus  measured  short-wave 
radiation  for  both  clear  and  cloudy  skies.  In  the  present  paper  the  amount 
of  reflected  short-wave  radiation  Qrr  was  estimated  by  fitting  a  polynomial 
to  average  values  between  Koberg's  clear  and  cloudy  sky  curves: 

QRR  -  °-i08QRI  •  ‘.IM-IO-'QJ,,  • 


(15) 
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Koberg  estimate*  that  thia  should  b«  accurate  to  within  10%.  The  amount  of 
energy  absorbed  by  the  water  is  then  calculated  by 


Under  condition*  where  Tw>Ta,  and  windspeed  is  low,  fog  generally 
occurs.  Such  conditions  are  more  likely  the  greater  the  difference  (Tw  - 
T4>  becomes.  Qp,  therefore,  is  commonly  appreciably  less  than  QCL  and 
contributes  only  slightly  to  the  total  energy  budget. 


The  net  heat  lost  from  the  water  surface  as  long-wave  radiation  is 
given  by 

*aB  *  Q.  *  Qar  *  Qb.  <17> 

where  Q*  is  the  incoming  long-wave  atmospheric  radiation,  Qar  is  the  re¬ 
flected  long-wave  atmospheric  radiation,  and  Qj,,  is  the  back  radiation 
from  the  water  surface.  The  Lake  Hefner  studies  (Anderson,  1954,  p.  96- 
9ft)  have  established  that  the  emissivity  of  water  is  0.  970  and  is  independent 
of  the  temperature  and  concentration  of  dissolved  solids.  Thus 

Qb8  *  0.970  eT^,  (18) 

where  <r  is  the  Stefan- Boltzmann  constant  (1.  171  x  10*Tcal  cm'1  day'1  ‘K"1) 
and  Tw  is  water  temperature  in  *K. 

Where  Qa  is  not  measured,  it  can  be  related  empirically  to  air  temper¬ 
ature,  vapor  pressure,  and  cloud  conditions  (Anderson,  1954,  p,  88-98), 
For  clear  skies  the  relation  initially  suggested  by  Brunt  is  used  with  the 
empirical  constants  determined  by  Anderson: 

Qa  «  (0.68  +  0.036  VeT)  eT«  (19) 

where  Ta  is  the  air  temperature  in  *K,  ea  is  the  vapor  pressure  of  the  air 
in  mb,  and  Qa  is  expressed  in  cal  cm"1  day"1 .  For  cloudy  skies  Qa  was 
estimated  by 


Qa  =  (4  +  bea>  *t; 


where 


a  *  0.740  +  0.025  C  exp  [( -1. 92)(10"4)H] 
b  =  (4.9)(10"J)  -  (5.  4)(10"«)  C  exp  ((-1.97){10-*)H] 


500  <  H  <  d. 
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Here  C  la  the  cloudinesa  expressed  In  tenths  and  H  la  the  cloud  height  in 
meters.  For  cloud  heights  less  than  S00  m,  the  cloud  height  is  taken  as 
a  constant  at  500  m  (Anderson,  1954).  Because  the  emissivity  of  water  is 
0.970,  its  long- wave  reflectivity  ia  0. 01.  so  that  Q4r  ■  0.  0]  Qa, 

If  weather  observations  are  available,  the  value  of  e#  can  be  computed 
from 


U'.a 

“nrcr 


(23) 


where  U  is  the  relative  humidity  and  eaa  is  the  saturation  vapor  pressure 
of  the  air  at  temperature  Ta.  The  value  of  ef».  which  is  a  function  of 
temperature  alone,  can  be  calculated  with  the  uoff-Gratch  formula  (List, 
1963.  p.  350). 

Qe  *  °H 

In  studies  of  this  type,  Q£  +  Of]  is  generally  estimated  by  using  an 
empirical  formula  to  calculate  either  Qg  or  Qjj,  and  then  using  the  Bowen 
ratio 


R 


to  find  the  other.  The  value  of  R  is  calculated  from 


(24) 


R 


■  (6. 1)(10“4) 


(25) 


where  P  is  atmospheric  pressure  (mb),  eaw  is  saturation  vapor  pressure 
at  the  temperature  of  the  water  (mb),  e*  is  vapor  pressure  of  the  air  (mb), 
and  other  symbols  are  as  previously  defined. 

Equation  25  is  derived  on  the  assumption  that  the  eddy  diffuaivities  of 
water  vapor  (i^)  and  heat  (K^)  ln  11  ir  ar e  equal.  Examination  of  p r evi ou s 
literature  reveals  uncertainty  about  the  conditions  under  which  this  assump¬ 
tion  holds  true.  Pasquill  (1949).  Veihmeyer  (1964)  and  Fritschen  (1965)  alt 
reported  satisfactory  results  using  eq  25  when  atmospheric  conditions  are 
neutral  or  stable.  Anderson  (1954,  p.  117)  also  concluded  as  a  result  of 
the  Lake  Hefner  studies  that  little  error  appears  to  be  introduced  in  the 
estimates  as  the  result  of  neglecting  the  effects  of  atmospheric  stability. 
Brutsaert  (1965,  p.  13)  cites  one  report  in  which  K  ■  i%/Kwwas  found  equal 
to  unity  regardless  of  atmospheric  stability,  and  another  where  the  ratio 
was  found  equal  to  unity  except  under  highly  stable  conditions,  where  it 
increased  to  2  or  3.  In  unstable  conditions  (TW»T4),  Pasquill  (1949)  re¬ 
ported  that  K  could  be  as  large  as  3.  A  similar  increase  in  the  value  of  K 
with  an  increase  in  (Tw  -  T.)  was  noted  by  Roll  (1965,  Fig.  61).  However, 
Rimsha  and  Donchenko  (1957),  whose  work  is  discussed  more  fully  below, 
successfully  used  the  Bowen  ratio  with  K  ■  1  to  calculate  wintertime  heat 
losses  from  calorimeters  and  open  river  reaches,  with  (Tw  -  Ta)  as  great 
as  30  C.  In  view  of  these  latter  studies,  it  seems  reasonable  to  use  the 
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Bowen  ratio  in  calculating  Qg  and  Qu  under  the  conditions  of  interest  here. 
Since  it  enters  the  theoretical  derivation  of  R  as  a  product  term,  the  lack 
of  information  as  to  its  true  vaiije  can  introduce  considerable  uncertainty 
in  the  estimates  of  or  Qg,  and  hence  in  the  total  heat  loss  rate.  The 
possibilities  of  such  error  must  be  kept  in  mind  in  situations  such  as  dis- 
cussed  herein. 

Using  the  Bowen  ratio,  two  approaches  to  calculating  Qg  and  Qh  were 
attempted.  The  first  involves  using  a  simple,  Well-accepted  formula  for 
calculating  daily  evaporation  rates  (Kohler,  1954). 

E  «  f  (0.  U25){10‘*  )  4  <1.007X10-*  )  vj  <##w  -  ea>  (26) 

where  E  is  evaporation  rate  (cm  day*1),  v4  is  average  wind  velocity  at  the 
4-m  level  (m  sec"1 ),  efrW  ia  saturation  vapor  pressure  at  the  temperature 
of  the  water  surface  (mb),  and  e4  is  the  vapor  pressure  of  the  air  at  the 
Z-m  level  (mb).  Multiplication  of  eq  26  by  suitable  constants  allows  the 
heat  lost  as  a  result  of  evaporation  to  be  calculated  by 

Qe  «  <3.  135  +  6.011  vj(e#w  -  e^)  (27) 

where  Qg;  is  in  cal  cm'1  day”1 ,  va  ia  in  m  sec“*  ,  and  esw  and  ca  are  in  mb. 
Assuming  that  eq  27  is  applicable  within  small  limits  of  error  when  2-m 
wind  velocities  are  used,  and  that  P  =  1000  mb,  its  multiplication  by  eq  2  5 
yields'* 

Qh  =  (1.91  4  3.67vm)(Tw  -  Ta).  (28) 

The  second  approach  uses  the  equations  developed  by  Rimsha  and 
Donchenko  (1957),  who  assumed  that  Q^j  could  be  calculated  by  a  relation 
of  the  form  of  eq  28,  and  determined  the  wind  function  parameters  by 
least  squares  analyses  of  data  from  calorimeters  exposed  adjacent  to 
natural  ice-free  reaches.  The  principal  difference  in  the  two  approaches 
lies  in  the  fact  that  the  "constant"  in  the  wind  function  was  found  by  Rimsha 
and  Donchenko  to  vary  with  (Tw  -  Ta).  Thus,  after  adjustment  for  appli¬ 
cation  to  natural  reaches  (based  on  actual  comparison  of  heat  losses  from 
rivers  and  calorimeters),  they^yrote _ _ _ 

°H  =  (kn  +  3-9va><Tw  -  Ta)  <*9> 

where  va  is  2-m  windspeed  in  m  sec-1  and  k„  is  empirically  found  to  be 

kn  *  8.0  ♦  0.35(TW  -  Ta).  (30) 


*  If  changes  are  made  in  eq  26,  27  and  28  to  correct  for  the  substitution  of 
average  wind  velocities  determined  at  the  2-m  level  in  equations  based  on 
4-m  wind  velocities,  these  equations  change  to 


E  - 

[  (0.  525)(1  0-*)  .  (1 . 22 9)(1  0-*)  vaJ  (e,  w-ca) 

(26a) 

Qe 

(3.135  i  7.  3 3  va)(e,w  -ea) 

(27a) 

«h 

(1.91  4  4.47  va)(Tw  ’  Ta). 

(28a) 

and 
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0  10  _ » _ 30 _ 40  _ 

VT0  <*C> 

Figure  l.  (Qg  ♦  Qj-f)  vs  {Tw  “  Ta)  at  va  =  0,  5,  and  1 0  m  sec'1, 

U  s  lOOft.  The  different  curves  are  calculated  by  the  Russian 
winter  equation  and  by  the  Kohler  equation  method. 


The  parameter  kn  thus  represents  the  portion  of  heat  transfer  due  to  free 
convection.  The  Bowen  ratio,  assuming  H  *  iOOO  mb,  was  then  applied  by 
Rimsha  and  Donchenko  (1957)  to  eq  29  to  give 

QE  •  (» .  56  kn  *  6.  08  va)  (efW  *  *»)<  0*) 

it  is  interesting  to  note  that  the  wind  coefficient  in  cq  31  is  very  close  to 
that  found  by  Kohler  (1954,  eq  1 7),  Thu*  cq  27  is  a  special  case  of  eq  31, 
which  is,  according  to  eq  30,  appropriate  for  (Tw  *  Ta)  <  1C. 

Future  2  compares  Q|.j  *  Qr  as  calculated  by  eq  27  and  28  with  that 
calculated  by  eq  29  and  3t,  ana  illust rates  the  large  differences  which  art.  '* 
front  introduction  of  a  term  accounting  for  natural  convection.  The  dif¬ 
ferences  are,  as  expected,  most  pronounced  at  low  windspeeds. 

It  is  recognized  that  eq  £6-3  1  are  simplifications,  and  that  the  con¬ 
stants  in  each  actually  depend  upon  such  factors  as  wind  sheltering,  fetch 
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length*,  end  the  amplitude  and  (rtqutPty  of  wind  v^rwUuha  «t  lU  nt«»  tut 
which  they  were  determined.  Some  field  tests  are  described  below,  but 
many  more  are  needed  before  the  general  applicability  of  one  or  the  other 
art  of  equation*  ia  demonstrated. 

Figure  3  shows  the  relative  importanre  nf  Op,  Qg  and  under 
differing  meteorological  conditions  using  eq  29,  30  ai.tl  31, 


Qs 

The  heat  loat  aa  a  reault  of  tnow  falling  into  the  water  Qg  (cal  cm"4 
day*1)  can  be  estimated  by  using  the  relation 

Qs  *  A(L  +  Ct(Tw  -  Ta)l  (W) 

where  A  ia  the  snow  accumulation  rate  (g  cm*4  day*1 ),  L  la  the  latent  heat 
of  fusion  of  ice  (cal  g'1),  and  C{  is  the  heat  capacity  of  ice.  If  possible, 

A  should  be  taken  from  direct  measurements  of  anowfall  rate.  In  many 
cases  these  do  not  exist,  but  estimates  of  visibility  will  be  available,  and 
a  satisfactory  estimate  may  be  possible  using  a  relationship  such  as  that 
suggested  by  Melior  (196b): 

A  «  7,85V*4*”*  (33) 


where  V  is  visibility  in  km.  Equation  3  3  was  found  to  hold  reasonably  well 
over  the  range  1  <  V  <  10  under  calm  conditions.  A  feel  for  the  relative 
Importance  of  Qs""in  contributing  to  the  value  of  Q*  can  be  obtained  by 
comparing  Figure  -i  with  Figure  3. 


An  exact  value  for  this  term  requires  a  knowledge  of  the  geothermal 
gradient  in  the  area  of  concern  and  the  thermal  dlffusivity  of  the  stream 
bed  material.  Lacking  this,  values  measured  under  similar  conditions 
may  be  found  in  the  literature.  The  only  reports  found  were  by  Incc  and 
Ashe  (1964),  who  indicated  that,  in  December,  Qq  at  5  cal  cm"4  day*1  for 
•ome^SwedisirriverB^snd  by^timaha  and  Donchenko  (1957)  who  gave  Qq  ~ 

10  to  20  cal  cm*4  day'1  for  some  Russian  rivers.  Both  of  these  values 
can  undoubtedly  be  considered  high  inasmuch  as  a  representative  heat  flow- 
value  for  the  Canadian  Shield  is  roughly  0.  1  cal  cm'4  day*1  (Lee  and  Uyeda, 
1965).  Thus  in  most  cases  it  is  safe  to  consider  Qq  negligible. 


°GW 

Fairly  detailed  knowledge  of  the  river  is  required  to  evaluate  this  term. 
Stream  discharge  hydrographs  can  be  analysed  to  provide  a  rough  estimate 
of  gains  or  losses  attributable  to  ground  water.  If  it  proves  necessary  to 
consider  ground  water,  the  mean  annual  air  temperature  at  the  site  should 
provide  a  reasonable  estimate  of  the  temperature  of  the  shallow  ground 
water.  Since,  if  ground  water  flow-  is  significant,  assumption*  of  constant 
discharge  arc  violated,  Qqw  i»  considered  negligible. 
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Fi«ure  3.  Qb.  Qe*  C>H*  *nd  a*  v*  (Tw  *  Ta)  at  v  0,  9,  and 
10  m  scc**j  U  -  50^j  clear-*ky  condition*.  Qg  and  Qh  calculated 
by  the  Ruatian  winter  equation  and  Bowen  ratio. 


] 


1  u  ir.'  1.  <Js  vs  (  [  w  -  T*)  at  va 


0. 
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The  heat  added  due  to  fluid  friction  can  be  calculated  using  a  relation 
given  by  Ince  and  Ashe  (1964,  p.  8): 

Q  =  8.  64  x  104  (34) 

F  J  v.' 

where  D  is  river  discharge  (err.*  sec"1  ),  y  is  weight  density  of  water  (dy 
cm’1),  S  is  the  slope  of  the  water  surface,  w  is  the  width  of  the  river  (cm) 
and  J  is  the  thermal  equivalent  (4.  ’87  x  1  07  dy  cm  cal"1  ).  Under  ir.o,t 
conditions,  Q  ^  is  insignificant. 


CALCULATION  OF  X 

The*  distance  downstream  from  the  site  ot  thermal  pollution  to  the 
cross-section  at  which  the  water  temperature  has  decreased  to  a  given 
value  is  found  by  substituting  the  chosen  expressions  for  the  heat  budget 
components  (eq  14-34)  in  ?q  11  and  integrating  between  suitable  limits.  It 
is  clear  that  if  the  meteorological  parameters  Ta,  va,  ea,  Qr^,  V,  C,  H 
and  P  can  be  considered  as  known  constants  over  the  reach  ot  river  of 
interest.  Q*  is  a  function  of  water  temperature  only.  Unfortunately,  the 
function 

[  dT  /Q* 

J  w 

cannot  be  integrated  in  closed  form  and  numerical  integration  must  be  used. 


Figure  5.  Typical  downstream  temperature  profile  below  a 
thermal  pollution  source. 


A  computer  program  was  prepared  to  provide  a  means  of  evaluating 
eq  1 1  by  using  Simpson's  rule.  The  program  was  prepared  so  that  calcula¬ 
tions  could  be  performed  for  any  combination  of  values  of  the  hydrometric 
and  meteorological  parameters.  The  water  temperature  at  the  site  of 
thermal  pollution  (x  =  0)  wa3  calculated  using  eq  12.  Simpson's  rule  was 
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then  used  to  compute  the  value  of  the  Integral  and  the  cot  responding  value 
of  x  in  0.  1  K  steps  downstream  until  some  cutoff  temperature  was  reached. 
Tn  the  problems  considered  in  this  paper  the  ,'utoff  temperature  was  i7i.  1  X, 
the  ice  point. 

Figure  3  shows  the  form  of  a  typical  temperature  vs  distance  curve. 
The  program  may,  of  course,  also  be  used  to  calculate  temperature  pro¬ 
files  downstream  from  thermal  pollution  sites  under  meteorological  con¬ 
ditions  such  that  no  ice  will  form,  as  long  as  the  river  is  cooling  at  the 
surface. 


"EXPERIMENTAL"  CHECKS 

It  is  difficult  to  state  the  accuracy  with  which  the  approach  outlined  in 
this  paper  will  permit  the  calculation  of  the  position  of  the  ice  front.  Be- 
s.des  the  difficulty  in  accurately  predicting  Q*,  it  is  obvious  that  an  ice 
sneet  will  not  form  exactly  at  the  point  where  the  water  temperature 
reaches  0C.  Under  conditions  of  intense  tutbulent  mixing,  frazil  ice  may 
form  initially.  Where  vertical  mijdng  is  weak,  the  theoretical  position  of 
the  ice  sheet  will  in  general  also  lie  upstream  from  the  actual  ice  front. 

It  should  be  emphasized  that  this  paper  calculates  only  the  steady-state 
ice-free  reach.  The  calculations  shov'd,  therefore,  be  most  representative 
when  the  weather  conditions  and  Qp  have  oven  relatively  constant.  When 
either  the  air  temperature  or  Qp  are  decreasing  with  time  the  calculations 
should  also  be  representative,  because  the  position  of  the  ice  front  should 
rapidly  adjust  to  the  new  steady-state  vafuc.  However,  the  length  of  the 
calculated  ice-free  reach  may  he  in  considerable  error  when  an  appreciable 
cold  peri'  is  followed  by  a  warming  trend.  During  the  cold  period,  ice 
will  accumulate  on  the  river  below  the  position  of  the  ice  front.  When  the 
air  temperature  rises,  some  of  this  ice  must  melt  as  the  length  of  the  ice- 
free  reach  increases  to  the  new  steady-state  value.  The  time  required  for 
this  ice  to  melt  will  cause  a  time  lag  in  achieving  the  new  steady-state 
position  of  the  ice  front. 

Because  of  these  difficulties  and  because  many  of  the  assumptions 
involved  in  writing  eq  11  are  not  strictly  true  in  a  real  situation,  it  is  highly 
— desi  rable  to  be  a  ble  to  -compare  field  measu  reme  nts-of-i  ee  -free  reac  hes 
below  existing  thermal  pollution  sites  with  computed  values.  Unfortunately, 
careful  field  observations  of  this  phenomenon  are  extremely  limited.  How¬ 
ever,  it  was  possible  to  obtain  unpublished  data  collected  by  the  Corps  of 
Engineers,  St.  Paul  District,*  regarding  the  length  of  ict-free  reaches 
below  the  Northern  State  Power  Compan/’s  Riverside  and  Highbridge  steam- 
power  generating  plants.  These  plants  are  located  on  the  Mississippi  River 
in  the  vicinity  of  Minneap'  ''s-St.  Paul  and  discharge  thermal  pollution 
directly  into  the  river.  a  test  area  this  general  reach  of  the  Mississippi 
hau  several  disadvantages:  several  other  small  plants  discharge  thermal 
pollution  and  small  tow  boats  frequently  pass  up  and  down  the  river  when 
ice  conditions  are  not  too  severe.  Both  these  factors  undoubtedly  disturb 
the  natural  development  of  an  ice  sheet.  The  extent  to  which  they  enlarge 
the  ice-free  reaches  below  the  Riverside  and  Highbridge  plants  is  unknown. 
Further,  at  a  point  5  km  below  the  Riverside  plant  there  is  a  lock  and  dam 


*  M.  Nelson,  personal  communication. 
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Tftbl*  t.  Observed  tnd  calculated  vaiot*  fot  the  area*  and  teoftba  of  the  ire-free  reaches  below  thermal 
pollution  *tt«a  aloof  the  Mississippi  River.  Minneapolis  St  Paul.  Minnesota. 
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complex,  and  *  4  km  below  the  Highbridge  plant  a  significant  amount  of 
thermal  pollution  is  introduced  by  the  Minneapo)is-St.  Paul  Sanitary  Dis- 
trict  plant.  Therefore  neither  of  these  locations  represents  an  ideal  test 
site  (i.e.,  a  location  where  an  isolated  thermal  pollution  source  introduces 
heat  into  a  river  that  possesses  relatively  constant  hydrometric  param¬ 
eters). 

The  total  area  of  open  water  below  these  two  plants  was  determined 
from  aerial  photographs  on  several  different  dates  during  the  winters  of 
1965  and  1966.  Because  the  open  water  -  ice  boundary  was  usually  quite 
irregular,  the  area  of  open  water  was  converted  to  an  effective  length  by 
dividing  it  by  a  representative  width  of  the  river.  Below  the  Riverside 
plant  the  full  width  of  the  river  was  used.  However,  below  the  Highbridge 
plant  an  effective  width  of  one-half  the  actual  width  was  utilized  because 
field  observations  indicated  that  thermal  mixing  was  confined  to  the  left 
half  of  the  river.  The  meteorological  data  used  are  from  the  U.S.  Weather 
Bureau  Station  at  the  Minneapolis -Si.  Paul  Airport  (average  daily  values), 
and  the  daily  average  thermal  pollution  values  (Qp,  cal  sec*1)  for  the 
Riverside,  Highbridge  and  Minneapolis -St.  Paul  Sanitary  District  plants 
were  provided  by  the  St.  Paul  District,  U.S.  Army  Corps  of  Engineers. 

These  data  are  tabulated  in  Appendix  C. 

The  question  immediately  arises  as  to  the  time  period  over  which  the 
meteorological  conditions  should  be  averaged  or,  stated  another  way,  how 
quickly  the  length  of  the  open  wat<  r  reach  responds  to  changes  in  meteoro¬ 
logical  conditions.  While  no  study  was  undertaken  to  resolve  this  matter, 
Bilello's  (1967)  study  of  Seneca  Lake,  Michigan,  suggests  that  near-surface 
water  temperature  trends  follow  air  temperature  trends  with  a  lag  of  about 
1  day.  The  work  of  Ince  and  Ashe  (1964)  on  the  St.  Lawrence  River  suggests 
a  similar  response  time.  In  most  cases  considered  in  the  present  paper,  the 
length  of  time  required  for  the  water  to  tr.. verse  the  ice-free  reach  is  less 
than  o.-  only  slightly  greater  than  1  day.  Therefore,  the  average  meteoro¬ 
logical  data  from  the  date  prior  to  the  date  of  ice  observation  were  used  in 
calculating  the  length  of  the  ice-free  reach. 

Table  I  presents  the  observed  ice*free  areas  below  the  thermal  pollution 
sites  and  "observed"  ice-  free  len  ths  caiculatec  from  these  areas.  Ice-free 
lengths  calculated  using  the  Kohler  equation  method  and  the  Russian  winter 
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OBSERVED  LENGTH  of  ICE-FREE  REACH,  km 

Figure  6.  Observed  vs  calculated  ice-free 
reaches.  Measurements  below  the  Highbridge 
Plant  are  shown  as  open  circles  and  measure¬ 
ments  below  the  Riverside  Plant  as  dots. 


equation  are  also  listed.  In  the  calculations  regarding  the  Highbridge 
plant,  the  calculated  open  reach  is  longer  than  the  distance  (6.4  km)  to  the 
next  large  thermal  pollution  source,  the  Minneapolis-St.  Paul  Sanitary  Dis¬ 
trict  plant.  Therefore,  to  calculate  the  total  length  of  the  open  reach,  the 
calculated  water  temperature  6.4  km  below  the  Highbridge  plant  was  taken 
as  the  natural  water  temperature  Tnat  when  the  thermal  pollution  associated 
with  the  Minneapolis-St.  Paul  Sanitary  District  plant  is  introduced.  The 

final  length  of  the  open  reach  produced  by  both  plants  is  then  simply  the _ 

length  of  fhe  calculated  reach  below  the  Minneapolis-St.  Paul  Sanitary  Dis¬ 
trict  plant  plus  6.4  km. 

As  would  be  expected  from  examining  Figure  i,  the  lengths  calculated 
by  the  Russian  winter  equation  are  always  shorter  than  those  calculated  by 
using  the  Kohler  equation  method.  They  are  also  clearly  in  better  agree¬ 
ment  with  the  observed  lengths.  Figure  6  presents  a  plot  of  the  observed 
lengths  of  the  ice-free  reaches  against  the  lengths  calculated  using  the  Russian 
winter  equation.  The  agreement  is  very  encouraging.  Although  there  is  a 
slight  tendency  for  the  calculated  lengths  to  be  longer  than  observed,  no 
consistent  pattern  is  apparent. 
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APPLICATIONS  TO  THE  ST.  LAWRENCE  SEAWAY 

The  ship  traffic  on  the  St.  Lawrence  Seaway  is  completely  stopped  by 
ice  for  3|  to  5  months  each  winter.  Particular  difficulty  is  encountered  in 
the  vicinity  of  the  lock  areas  and  ship  canals.  If  these  areas  could  be  kept 
ice-free  it  would  be  possible  to  move  shipping  into  the  Great  Lakes  during 
most  of  the  winter.  This  would  be  a  great  convenience  and  would  save  trans¬ 
portation  costs  of  several  million  dollars  per  year.  One  conceivable  way  to 
accomplish  this  is  to  use  the  excess  heat  produced  by  nuclear  pow«r  reactors 
to  remove  ice  at  critical  locations.  The  production  of  electricity  by  nuclear 
power  is  now  becoming  commercially  competitive  with  conventional  forms  of 
power  production.  In  addition,  the  exact  positioning  of  power  reactors  is  not 
important  as  long  as  there  is  a  major  river  nearby  that  can  be  used  as  a 
source  of  coolant  and  the  reactor  can  be  coupled  with  the  general  power  grid 
system  currently  in  use  in  eastern  Canada  and  the  United  States.  If  it  can 
be  shown  that  the  thermal  pollution  from  a  reactor  is  sufficient  to  keep  a 
significant  length  of  a  shipping  canal  ice-free  during  the  winter,  a  strong 
case  can  be  made  for  positioning  future  reactors  along  the  St.  Lawrence 
River. 

Two  types  of  nuclear  reactors  are  considered  here:  those  presently 
available  with  a  thermal  pollution  output  of  2.  9  x  10*  cal  sec"1  (~  1.2  x  I  0* 
watts)  and  those  >n  the  design  stage  with  a  thermal  output  of  8.  4  x  10*  cal 
sec*1  ('  3,  5  x  10*  watts).  As  an  estimate  of  extremely  severe  winter  con¬ 
ditions  for  the  general  seaway  area,  a  value  of  Ta  =  -17C  was  used.  This 
is  based  on  the  mean  monthly  minimum  temperatures  from  Ottawa  {January, 
-16.  7C;  February,  *15. 6C;  31  years  of  record)  and  Quebec  (January,  *16. 1C; 
February,  -15.  6C;  32  years  of  record).  A  windspeed  of  5  m  sec"*  was  as¬ 
sumed  as  representative,  based  on  January  and  February  data  from  Cornwall, 
Ontario.  This  windspeed  is  undoubtedly  high  because  in  the  general  St,  Law¬ 
rence  area  low  temperatures  are  usually  associated  with  relatively  calm 
conditions  and  pronounced  radiational  cooling.  A  representative  average 
temperature  for  the  seaway  area  during  January  and  February -D>—l-l-G-»--- 

_ As-an  example  the  ship'cahal  around  the  Montreal  area  is  considered: 

the  South  Shore  Canal  with  the  St.  Lambert  and  Cote  St.  Catherine  Locks, 
in  the  fall  and  early  winter  the  water  temperatu res  in  this  area  are  1  to  2 
degrees  lower  than  upstream  in  the  vicinity  of  Cornwall;  therefore,  the 
South  Shore  Canal  is  the  first  location  on  the  Seaway  that  encounters  signifi¬ 
cant  ice  problems.  The  average  width  of  the  canal  is  1  50  m,  the  depth  is 
9  .n  and  the  water  flow  velocity  0.  5  m  sec-*.  The  len,  th  of  the  South  Shore 
Canal  from  Montreal  Harbor  to  Lake  St.  Louis  is  roughly  32  km.  As  can  be 
seen  from  Figure  7  the  lengths  of  the  Canal  that  can  be  kept  ice-free  with 
Ta  -  -17C,  va  -  5  rr>  sec'1  and  C  IO  are  17  (Russian  winter)  and  25  km 
(Kohler*8  equation)  using  a  conventional  reactor.  With  a  larger  reactor  the 
calculated  lengths  are  50  and  75  km  respectively.  If  similar  estimates  are 
made  by  using  the  Russian  winter  equation  and  the  average  winter  temperature 
of  -l  1C,  values  of  27  and  75  km  are  obtained  for  the  conventional  and  design 
reactors  respectively.  It  is  obvious  that  even  a  conventional  reactor  would 
considerably  lengthen  the  shipping  season,  while  the  larger  reactor  would 
effectively  keep  the  canal  open  ail  winter. 

These  figures  indicate  the  definite  feasibility  of  using  nuclear  reactor 
power  plant  cooling  effluent  to  keep  selected  portions  of  the  St.  Lawrence 
Seaway  ice-free  for  much  of  the  winter.  The  prevention  of  ice  formation  in 
the  vicinity  of  locks  would  be  especially  beneficial. 
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Figure  7.  Length  of  ice-free  reach  as  afunction  of  air  temperature  (*C) 
and  wind  speed  (m  sec'1):  South  Shore  Canal,  St.  Lawrence  Seaway. 

Solid  curves  are  based  on  the  Russian  winter  equation  and  the  dashed 
curves  on  the  Kohler  equation  method. 

CONCLUSIONS 

_ The  general  method  far  calculating  the  downstream  temperature  profile _ 

of  a  river  and  the  length  of  reach  which  can  be  kept  ice-free  by  artificial  heat 
addition  was  stated,  but  not  developed,  by  Pruden  et  al.  (1954).  They  and 
previous  workers  did  not  develop  a  universally  applicable  expression  for  the 
heat  loss  rate  Q*  which  could  be  used  with  daily  averages  of  the  meteorologic 
variables.  The  present  paper  derives  the  basic  temperature-distance  rela¬ 
tion  from  fundamental  considerations  and,  by  use  of  a  computer  program, 
allows  the  use  of  several  of  the  many  previously  developed  empirical  and 
semiempirical  formulas  for  the  heat  budget  terms  comprising  Q*.  The  com¬ 
puter  program  thus  makes  possible  more  universally  applicable  estimates  of 
temperature  profiles  and  lengths  of  ice-free  reaches  than  were  heretofore 
practicable.  It  has  at  least  three  important  additional  advantages:  1)  tem¬ 
perature  profiles  and  lengths  of  ice-free  reaches  can  be  quickly  calculated 
for  any  combination  of  heats  of  pollution  and  meteorologic  and  hydrometric 
variables;  Z)  downstream  changes  in  meteorologic  and  hydrometric  variables 
can  be  accounted  for,  if  desirable,  by  dividing  the  river  into  any  number  of 
sub-reaches  and  treating  the  water  temperature  at  the  outlet  of  one  sub-reach 
as  the  temperature  at  the  beginning  of  the  next  sub-reach  downstream;  and 
3)  the  computer  program  can  be  easily  modified  to  permit  the  use  of  differ¬ 
ent  expressions  for  the  heat  budget  items  or  the  addition  of  other  heat  budget 
terms  which  may  be  significant  in  a  particular  case. 

The  important  problems  that  must  be  overcome  to  improve  the  "accura¬ 
cy'7  of  the  present  method  of  calculation  are  the  limited  understanding  of  the 
evaporative  and  convective  heat  loss  processes  when  an  unstable  vertical 
density  distribution  exists  over  open  water,  and  the  lack  of  a  general  method 
for  evaluating  the  vertical  and  especially  the  horizontal  mixing  of  heat  in  a 


THE  EFFECTS  OF  THERMAL  POLLUTION  ON  RIVER  ICE  CONDITIONS  41 

river.  It  is,  however,  encouraging  that  observations  on  the  Mississippi  River 
at  Minneapolis-St.  Paul  show  reasonable  agreement  with  computed  lengths  ot 
ice-free  reaches.  Better  tests  of  the  method  are  highly  desirable. 

Calculations  based  on  representative  values  of  the  meteorologic  and  hydro¬ 
metric  variables  for  the  South  Shore  Canal  portion  of  the  St.  Lawrence  Seaway 
funnest  that  significant  reaches  of  that  waterway  could  be  kept  tee-free  bv  the 
proper  placement  of  nuclear  reactors.  Significant  lengthening  of ^  the  shipping 
season  with  concomitant  economic  advantages  would  thus  be  feasible. 

Although  the  current  computer  program  is  not  designed  <0  handle  situation* 
where  the  air  temperature  is  warmer  than  the  water  temperature,  it  can  easily 
be  modified  to  include  such  situations.  It  this  were  done,  the  program  would  be 
useful  in  treating  such  problems  as  the  cumulative  effect  of  the  positioning  ot  a 
number  of  nuclear  reactors  along  the  Columbia  River  on  the  temperature  of  the 
river  (Geraghty,  1967). 
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APPENDIX  A:  COMPUTER  PROGRAM  FOR  CALCULATING  WATER  23 
TEMPERATURES  BELOW  THERMAL  POLLUTION  SITES 


Symbology 

The  postscript  ---IN  denotes  the  Initial  value  of  the  variable  •**,  the 
postscript  *--DEL  indicates  the  incremental  value  of  — ,  and  the  postscript 
— F  indicates  the  final  value.  The  variable  — INC  denotes  the  initial  value 
of  —  In  degrees  centigrads,  ate.  The  units  to  be  used  for  the  variables 
are  as  listed  in  the  Liet  of  Symbols.  Because  the  *-*DEL  values  are  always 
positive,  — IN  should  bs  a  smaller  value  than  — F  (e.g.,  VAIN  =  #. ,  VADEL 
=  S.#,  VAF  -  2#.#.,  or  TAINC  «  ‘5#,#,  TADEL  *  If.#,  TAFC  * 
which  correspond  to  TAIN  =  223.2,  TADEL  -  1#.#,  TAF  =  263.2).  1NT  is 
the  number  of  intervale  used  in  the  Simpson's  rule  calculations  a  value  of 
INT  *1#  was  found  satisfactory. 

Senas  Switch  Instructions 

Sense  Switch  1 

OFF:  POLL  *  POLL  ♦  POLLDEL  (FOLLDEL  is  uted  and  there 
is  a  constant  additive  increment). 

ON:  POLL  =  1#.  *  POLL  {the  program  increment#  POLL  by 
powers  of  !#).  In  this  case  any  fixed  point  number  ie>  put  in  the  POLLDEL 
site  on  the  data  tape. 

Sense  Switch  2 

OFF:  Not  snowing  (this  option  is  given  beca  ise  decreased 
visibility  does  not  necessarily  indicate  snow). 

ON:  Snowing. 

Sense  Switch  3 

OFF:  The  Russian  equations  (eq  £9-31)  are  used  to  calculate 
QE  and  QH. 

ON:  The  Kohler  equations  (oq  27-28)  are  used  to  calculate 
QE  and  QH. 

Data  Input 

The  data  are  introduced  in  the  order  and  format  specified  in  the  program. 
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APPENDIX  B:  COMPUTER  PROGRAM  FOR  CALCULATING  HEAT 
BUDGET  COMPONENTS  AND  TOTAL  HEAT  LOSS 
RATES  FROM  A  COOLING  RIVER 


Symbology 

Initial,  incremental,  and  final  valuea  of  variables  are  designated  as 
described  in  Appendix  A. 

Sense  Switch  Instructions 

Sense  switches  are  used  as  in  Appendix  A,  with  the  exception  that  Sense 
Switch  1  is  not  used. 

Data  Layout 

The  data  are  introduced  in  the  order  and  format  specified  in  the  pro* 
gram.  The  FUNCTION  VAP(Z)  subroutine  used  in  this  program  is  identical 
to  that  used  in  Appendix  A. 
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An  attempt  is  made  to  calculate  the  length  of  the  ice -free  reach  which  develops 

during  the  winter  below  a  thermal  pollution  site  on  a  river. _ A  differential  equation 

for  the  steady  state  heat  balance  of  a  volume  element  of  a  river  is  developed,  which 
leads  to  the  expression 


where  x  is  distance  downstream  from  the  pollution  site  to  the  cross  section  where 
the  water  temperature  equals  Twx,  Twq  is  water  temperature  at  x  equals  zero,  Q* 

is  rate  of  heat  loss  from  the  water  surface,  and  Cx  is  a  constant  which  includes  flow 
velocity  and  depth.  The  value  of  x  at  Twx  equals  0*C  is  taken  as  the  length  of  the 
ice-free  reach.  Q*  is  the  sum  of  heat  losses  due  to  evaporation,  convection,  long- 
and  short-wave  radiation,  and  other  processes,  each  of  which  is  evaluated  by  an 
empirical  or  theoretical  expression.  The  two  principal  limitations  in  accurately 
calculating  downstream  temperature  changes  are  related  to  difficulties  in  evaluating 
the  degree  of  lateral  mixing  in  natural  rivers  and  the  convective  and  evaporative  heat 
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losses  under  unstable  atmospheric  conditions. 
Observations  of  lengths  of  ice-frce  reaches  on 
the  Mississippi  River  are  in  good  agreement  with 
the  calculated  values.  Significant  portions  of  the 
St.  Lawrence  Seaway  can  be  kept  ice-free  by  the 
installation  of  nuclear  reactors  at  appropriate  lo¬ 
cations.  i 


